//
//  $Author:: chrisandsusan3587@gmail.com                  $
//    $Date:: 2012-07-27 02:57:00 +0000 (Fri, 27 Jul 2012) $
//     $Rev:: 411                                          $
//

//  CPE 790
//  Chris Davis
//  Final Project

//  averagekernel.cu
//  Kernel file that will average blocks of data.

#include "project.h"


__global__ void GPUSum(unsigned int * inputData, unsigned int * outputData, int NumberOfElementsToAverage)
{
	unsigned int partialsum = 0;
	unsigned int temp = blockDim.x;

	int idx = blockIdx.x * (NumberOfElementsToAverage / blockDim.x) + (threadIdx.x * ((NumberOfElementsToAverage / blockDim.x)/temp));

	for ( int i = idx; i < idx + ((NumberOfElementsToAverage / blockDim.x)/temp); i++)
	{
		partialsum += inputData[i];
	}

	outputData[(blockIdx.x * blockDim.x) + threadIdx.x] = partialsum;
	//__syncthreads();
	//atomicAdd(&outputData[0], partialsum);
}

double computeGPUAverage(unsigned int * cube,
                      unsigned long Xsize,
                      unsigned long Ysize,
                      int StartingWavelength,
                      int WavelengthsToAverage,
                      int blockSize,
                      int gridSize)
{
	int NumberOfElementsToAverage = 0;
	int WavelengthOffset = 0;
	double averagedResult = 0.0;
	unsigned int * d_deviceCube;
	unsigned int * d_resultCube;
	unsigned int * h_resultCube;
	//int blockSize;
	//int gridSize;
	int totalSum = 0;

	NumberOfElementsToAverage = (Xsize*Ysize) * WavelengthsToAverage;
	WavelengthOffset = (Xsize*Ysize) * StartingWavelength;

	//blockSize = Xsize * 2;
	//gridSize  = Xsize * 2;

	h_resultCube = (unsigned int *) malloc ( gridSize * blockSize * sizeof(unsigned int) );

	cudaMalloc( (void **)&d_deviceCube, NumberOfElementsToAverage * sizeof(unsigned int) );
	cudaMalloc( (void **)&d_resultCube, gridSize * blockSize * sizeof(unsigned int) );

	cudaMemcpy( d_deviceCube, &cube[WavelengthOffset], NumberOfElementsToAverage*sizeof(unsigned int), cudaMemcpyHostToDevice);
	cudaMemcpy( d_resultCube, h_resultCube, NumberOfElementsToAverage*sizeof(unsigned int), cudaMemcpyHostToDevice);

	GPUSum<<<gridSize, blockSize>>>(d_deviceCube, d_resultCube, NumberOfElementsToAverage);

	cudaMemcpy( h_resultCube, d_resultCube, gridSize * blockSize * sizeof(unsigned int), cudaMemcpyDeviceToHost );

	for(int i = 0; i < (gridSize * blockSize); i++)
	{
		totalSum += h_resultCube[i];
	}
	averagedResult = (double)totalSum / (double)NumberOfElementsToAverage;

	cudaFree(d_deviceCube);
	cudaFree(d_resultCube);

	free(h_resultCube);

	return averagedResult;
}